Load in data

library(tidyverse)
library(readxl)
library(ggplot2)
library(choroplethr)
library(choroplethrMaps)
library(gridExtra)
library(extracat)
library(viridis)
# Read in data
fea <- read_csv('fea_04062017.csv')

Analysis of Data Quality

# Missing data for Access and SNAP variables
missing_main <- select(fea, FIPS, State, County, LACCESS_POP10, PCT_LACCESS_POP10, LACCESS_LOWI10,  PCT_LACCESS_LOWI10,LACCESS_CHILD10, PCT_LACCESS_CHILD10, LACCESS_SENIORS10, PCT_LACCESS_SENIORS10, LACCESS_HHNV10, PCT_LACCESS_HHNV10, REDEMP_SNAPS08,  REDEMP_SNAPS12, PCH_REDEMP_SNAPS_08_12, PCT_SNAP09, PCT_SNAP14, PCH_SNAP_09_14, PC_SNAPBEN08,   PC_SNAPBEN10,   PCH_PC_SNAPBEN_08_10)
visna(missing_main)

# Missing data for Food Insecurity variables
missing_main_2 <- select(fea, State, County, FOODINSEC_00_02, FOODINSEC_07_09, FOODINSEC_10_12, CH_FOODINSEC_02_12, CH_FOODINSEC_09_12, VLFOODSEC_00_02, VLFOODSEC_07_09, VLFOODSEC_10_12, CH_VLFOODSEC_02_12   , CH_VLFOODSEC_09_12, FOODINSEC_CHILD_01_07, FOODINSEC_CHILD_03_11)
visna(missing_main_2)

# Missing data for Prices, Taxes, and Restaurant variables
missing_main_3 <- select(fea, State, County, MILK_PRICE10, FOOD_TAX11,  FFRPTH07,   FFRPTH12,   PCH_FFRPTH_07_12,   FSRPTH07,   FSRPTH12,   PCH_FSRPTH_07_12,   PC_FFRSALES02, PC_FFRSALES07, PC_FSRSALES02)
visna(missing_main_3)

# Missing Data for Income and Poverty Variables
missing_main_4 <- select(fea, State, County, MEDHHINC10, POVRATE10, CHILDPOVRATE10)
visna(missing_main_4)

Main Analysis

# Quick data cleaning
# Remove PR - only look at 50 states + DC
fea <- filter(fea, State != "PR")
# For county plotting purposes
fea$region <- as.numeric(fea$FIPS)

As our project focuses on Food Insecurity in the US, the first variable(s) we will begin to analyze are those related to food insecurity. The documentation defines household food insecurity as those households which “were unable, at times during the year, to provide adequate food for one or more household members because the household lacked money and other resources for food.” The documentation adds that this inadequacy was in both “quality and variety of foods”. Households with very low food insecurity are categorized as households where “food intake of one or more members was reduced and eating patterns were disrupted at times during the year because of insufficient money and other resources for food.”

Note: Food insecurity data is at the state level

fea %>% 
  select(-region) %>% 
  group_by(State) %>% 
  summarise_all(mean) -> 
  state_data

# Bring in state abbreviation to name mapping
data(state)
state_abbrevs <- tibble(abb = state.abb, region = tolower(state.name))
state_abbrevs <- add_row(state_abbrevs, abb = 'DC', region = 'district of columbia')

state_data <- left_join(state_data, state_abbrevs, by = c('State' = 'abb'))
# Set manual scale for both graphs
# Remove state abbreviation labels
breaks = seq(6, 21, by = 3)

c = StateChoropleth$new(mutate(state_data, value = cut(state_data$FOODINSEC_07_09, breaks = breaks)))
c$title = "Household Food Insecurity (%), 2007-2009"
c$legend = "Food Insecurity %"
c$set_num_colors(length(breaks) - 1)
c$set_zoom(NULL)
c$show_labels = FALSE
c$ggplot_polygon = geom_polygon(aes(fill = value), color = 'gray20')
state_insec_0709 = c$render()

c = StateChoropleth$new(mutate(state_data, value = cut(state_data$FOODINSEC_10_12, breaks = breaks)))
c$title = "Household Food Insecurity (%), 2010-2012"
c$legend = "Food Insecurity %"
c$set_num_colors(length(breaks) - 1)
c$set_zoom(NULL)
c$show_labels = FALSE
c$ggplot_polygon = geom_polygon(aes(fill = value), color = 'gray20')
state_insec_1012 = c$render()

state_insec_0709

state_insec_1012

#par(mfrow = c(1, 2))
#state_choropleth(mutate(state_data, value = cut(state_data$FOODINSEC_07_09, breaks = breaks)), num_colors = length(breaks) - 1, legend = 'Food Insecurity %', title = "Household Food Insecurity (%), 2007-2009")
#state_choropleth(mutate(state_data, value = cut(FOODINSEC_10_12, breaks = breaks)), num_colors = length(breaks) - 1, legend = 'Food Insecurity %', title = "Household Food Insecurity (%), 2010-2012")

The below charts explore average food insecurity for 2007-2009 and 2010-2012 by state. Here we see the range of percent of households experiencing food insecurity is increasing: from 07-09 the range was 6.7-17.7%, whereas from 10-12 the range was 8.7-20.9%.

# Cleveland dot plot theme
theme_dotplot1 <- theme_bw(10) +
    theme(axis.text.y = element_text(size = rel(.8)),
          axis.ticks.y = element_blank(),
          axis.title.x = element_text(),
          panel.grid.major.x = element_blank(),
          panel.grid.major.y = element_line(size = 0.5),
          panel.grid.minor.x = element_blank())
insec_bar <- state_data %>% select(State, FOODINSEC_07_09, FOODINSEC_10_12)

ggplot(insec_bar, aes(x = reorder(State, FOODINSEC_07_09 / 100), y = FOODINSEC_07_09 / 100)) + 
  scale_y_continuous(labels = scales::percent) +
  geom_point() + 
  theme_dotplot1 + 
  coord_flip() +
  ggtitle("Household food insecurity by state (3-year avg) - 2007-2009") +
  ylab("Household food insecurity (3-year avg)") + xlab("")

ggplot(insec_bar, aes(x = reorder(State, FOODINSEC_10_12 / 100), y = FOODINSEC_10_12 / 100)) + 
  scale_y_continuous(labels = scales::percent) + 
  geom_point() + 
  theme_dotplot1 + 
  coord_flip() +
  ggtitle("Household food insecurity by state (3-year avg) - 2010-2012") +
  ylab("Household food insecurity (3-year avg)") + xlab("")

##### Change colors?
insec_bar_tidy <- insec_bar %>% gather("Year", "Value", -State)

ggplot(insec_bar_tidy, aes(x = reorder(State, Value), y = Value / 100, color = Year)) + 
  scale_y_continuous(labels = scales::percent) + 
  geom_point() + 
  theme_dotplot1 + 
  coord_flip() +
  ggtitle("Household food insecurity by state (3-year avg) - 2007-09 vs. 2010-12") +
  ylab("Household food insecurity (3-year avg)") + xlab("") + 
  scale_color_manual(labels = c('2007-09', '2010-12'), values = c('#F8766D', '#00BFC4'))

##### Change colors?

# Density plot of food insecurity, comparing years
med_0709 <- median(filter(insec_bar_tidy, Year == 'FOODINSEC_07_09')$Value / 100)
med_1012 <- median(filter(insec_bar_tidy, Year == 'FOODINSEC_10_12')$Value / 100)

ggplot(insec_bar_tidy) + 
  geom_density(aes(Value / 100, fill = Year, color = Year), alpha = 0.4, adjust = 0.7) +
  scale_x_continuous(labels = scales::percent) + 
  xlab('Food Insecurity') + 
  ggtitle('Household food insecurity (3-year avg), with Medians') + 
  scale_fill_manual(name = 'Year', labels = c('2007-09', '2010-12'), values = c('#F8766D', '#00BFC4')) + 
  scale_color_manual(name = 'Year', labels = c('2007-09', '2010-12'), values = c('#F8766D', '#00BFC4')) + 
  geom_vline(xintercept = med_0709, color = '#F8766D', alpha = 0.8) + 
  geom_vline(xintercept = med_1012, color = '#00BFC4', alpha = 0.8) + 
  geom_text(aes(x = med_0709, y = -0.6, label = paste(med_0709*100, '%', sep = ''), hjust = 1.1), color = '#f86d6d') + 
  geom_text(aes(x = med_1012, y = -0.6, label = paste(med_1012*100, '%', sep = ''), hjust = -0.1), color = '#0091c2')

We’ll use this information to create subsets of the US, with Mississippi, Arkansas, Texas, and Alabama comprising the “high” food insecurity states, and North Dakota, Virginia, New Hampshire, and Minnesota making up the “low” food insecurity states.

high_insecurity_states <- c('MS', 'AR', 'TX', 'AL')
low_insecurity_states <- c('ND', 'MN', 'WI', 'MA')
fea_high <- filter(fea, State %in% high_insecurity_states)
fea_high$State <- paste('high', sep = '_', fea_high$State)
fea_low <- filter(fea, State %in% low_insecurity_states)
fea_low$State <- paste('low', sep = '_', fea_low$State)

Let’s explore variables that may be associated with food insecurity

Median household income

We expect median household incomes to be low in areas of high food insecurity and high in areas of low food insecurity. The below graphs look at county median household income across the US, and focuses in on states with the highest and lowest food insecurity rates.

continent_states <- unique(state_data$region)
continent_states <- continent_states[continent_states != 'alaska' & continent_states != 'hawaii']
blue_color_scale <- c('#084594', '#2171b5', '#4292c6', '#6baed6', '#9ecae1', '#c6dbef', '#eff3ff')
breaks = c(20, 33, 36, 39, 42, 45, 51, 120)
county_choropleth(mutate(fea, value = cut(MEDHHINC10/1000, breaks = breaks)), title = 'Median Household Income by County (2010)', state_zoom = continent_states) + 
  scale_fill_manual(values = blue_color_scale, name = 'Income ($1000s)')

#county_choropleth(mutate(fea, value = MEDHHINC10), legend = 'Income Quantile', title = 'Median Household Income by County (2010)', state_zoom = continent_states)
ggplot() + 
  geom_density(data = fea_high, aes(MEDHHINC10, color = State, fill = State), alpha = 0.1) + 
  geom_density(data = fea, aes(MEDHHINC10, color = 'US', fill = 'US'), alpha = 0.1) + 
  xlab('Median Household Income (County)') +
  ggtitle('Median Household Income (County): \n Highest Food Insecurity States & US (2010)') + 
  scale_x_continuous(labels = scales::dollar)

ggplot() + 
  geom_density(data = fea_low, aes(MEDHHINC10, color = State, fill = State), alpha = 0.1) + 
  geom_density(data = fea, aes(MEDHHINC10, color = 'US', fill = 'US'), alpha = 0.1) + 
  xlab('Median Household Income (County)') +
  ggtitle('Median Household Income (County): \n Lowest Food Insecurity States & US (2010)') + 
  scale_x_continuous(labels = scales::dollar)

As expected, county median household incomes for highly food insecure states are lower than for the US as a whole, while low food insecurity states have incomes distributed higher than for the US in general.

# Median HH income vs food insecurity
ggplot(fea, aes(x = MEDHHINC10, y = FOODINSEC_10_12 / 100)) + geom_point(alpha = 0.4) + 
  scale_x_continuous(labels = scales::dollar) + 
  scale_y_continuous(labels = scales::percent) + 
  xlab('Median Household Income') + 
  ylab('Average Food Insecurity') +
  ggtitle('State Average Household Food Insecurity (2010-12) \n vs County Median Household Income (2010), with Median of Counties') + 
  geom_vline(xintercept = median(fea$MEDHHINC10, na.rm = TRUE), linetype = 'longdash')

Poverty

Similarly, we can examine the relationship between the county-level poverty rate and state-level food insecurity rate. The poverty rate is defined as the percentage of county residents with household income below the poverty threshold, and the child poverty rate is defined as the percentage of county residents under age 18 living in households with income below the poverty threshold.

# We can use poverty and median hh income as proxies for each other - strong inverse relationship
ggplot(fea, aes(x = MEDHHINC10, y = POVRATE10 / 100)) + geom_point(alpha = 0.4) + 
  scale_x_continuous(labels = scales::dollar) + 
  scale_y_continuous(labels = scales::percent) + 
  xlab('Median Household Income') + 
  ylab('Poverty Rate') + 
  ggtitle('Poverty Rate vs Median Household Income (2010), by County')

Plotting state-level against county-level data makes understanding the true relationship between food insecurity and poverty difficult, especially because variance of poverty rates can be high within a state. Nevertheless, we can see that there is somewhat of a positive association between food insecurity with poverty and child poverty rates. Also, the states with highest food insecurity tend to have the majority of county poverty rates above the US county median.

Access to food

Part of our initial hypothesis was that food insecurity is related to low access to food. The dataset defines low access to store (%) as the “percentage of people in a county living more than 1 mile from a supermarket, supercenter or large grocery store if in an urban area, or more than 10 miles from a supermarket or large grocery store if in a rural area.” Is low store access more of an issue in counties with high food insecurity?

ggplot() + 
  geom_boxplot(data = fea_high, aes(x = State, y = PCT_LACCESS_POP10 / 100)) + 
  geom_boxplot(data = fea_low, aes(x = State, y = PCT_LACCESS_POP10 / 100)) + 
  geom_boxplot(data = fea, aes(x = 'US', y = PCT_LACCESS_POP10 / 100)) + 
  scale_y_continuous(labels = scales::percent) + 
  ylab('Low Access to Store (County)') + 
  ggtitle('Low Access to Store (County), with US Median (2010)') + 
  geom_hline(yintercept = median(fea$PCT_LACCESS_POP10 / 100, na.rm = TRUE), linetype = 'longdash') +
  geom_rect(aes(xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf), alpha = 0.1, fill = "red") + 
  geom_rect(aes(xmin = 5.5, xmax = 9.5, ymin = -Inf, ymax = Inf), alpha = 0.2, fill = "cyan3") + 
  xlim(c(unique(fea_high$State), 'US', unique(fea_low$State))) + 
  theme_bw()

From the above graph with multiple box plots, low access to store does not appear to be associated with food insecurity. In fact, at least half of the counties in AL, AR, and MS have low store access numbers lower than the US county median, while a higher proportion of people in MN and ND, where food insecurity is low, live far from stores. This suggests that physical access to stores is not a main factor related to food insecurity.

Other

# Define buckets of sales tax
breaks <- c(-1, 0, 2, 4, 6, 8)
state_data %>% 
  select(region, FOOD_TAX11) %>% 
  mutate(value = cut(FOOD_TAX11, breaks = breaks)) -> 
  plot_df
# Change first level to 0
levels(plot_df$value)[levels(plot_df$value) == '(-1,0]'] <- '0'

state_choropleth(plot_df, legend = 'Sales Tax %', title = "General Food Sales Tax (2011)")

It’s ironic that 3 of the 4 states (AL, AR, MS) with the highest food insecurity rates levy a general food sales tax, while most states in the US do not have a general food sales tax. Is this tax related to food insecurity, and if so, is this a chicken-and-egg problem?

ggplot() + 
  geom_point(data = fea_high, aes(x = State, y = MILK_PRICE10)) + 
  geom_point(data = fea_low, aes(x = State, y = MILK_PRICE10)) + 
  geom_boxplot(data = fea, aes(x = 'US', y = MILK_PRICE10)) + 
  ylab('Price of low-fat milk / national average') + 
  ggtitle('Regional price of low-fat milk / national average (2010)') + 
  scale_x_discrete(limits = sort(c('US', unique(fea_low$State), unique(fea_high$State)), decreasing = T)) +
  geom_rect(aes(xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf), alpha = 0.1, fill = "red") + 
  geom_rect(aes(xmin = 5.5, xmax = 9.5, ymin = -Inf, ymax = Inf), alpha = 0.2, fill = "cyan3") + 
  geom_hline(yintercept = 1, alpha = 0.3) + 
  xlim(c(unique(fea_high$State), 'US', unique(fea_low$State))) + 
  theme_bw()

##### Change colors?

# Farm analysis
store_farm_09 <- fea %>% select(FIPS, State, County, GROCPTH07, SUPERCPTH07, FMRKTPTH09, FOODINSEC_07_09)
store_farm_09_tidy <- store_farm_09 %>% gather("StoreType", "ValuePer1000", -FIPS, -State, -County,-FOODINSEC_07_09)

##By food insecurity level, store types (07-09)
ggplot(store_farm_09_tidy, aes(x = FOODINSEC_07_09 / 100, y = ValuePer1000, color = StoreType, fill = StoreType)) + 
  geom_col() + 
  scale_x_continuous(labels = scales::percent) + 
  ggtitle('Number of Stores per 1,000 pop') + 
  xlab('State Food Insecurity Rate (2007-09)') + 
  ylab('Stores per 1,000 pop') + 
  scale_color_manual(name = 'Store Type', labels = c('Farmers\' markets (2009)', 'Grocery stores (2007)', 'Supercenters & \n club stores (2007)'), values = c('#F8766D', '#00BA38', '#619CFF')) + 
  scale_fill_manual(name = 'Store Type', labels = c('Farmers\' markets (2009)', 'Grocery stores (2007)', 'Supercenters & \n club stores (2007)'), values = c('#F8766D', '#00BA38', '#619CFF'))

Diving into states

Here we dive into the four states with the highest and lowest food insecurity, according to the 2010-2012 average rate. The states with the highest food insecurity are: MS, AR, TX, AL; and the 4 states with the lowest food insecurity are: ND, VA, NH, MN. Let’s see what differentiates these high and low food insecurity states.

high_low_states <- c('mississippi', 'arkansas', 'texas', 'alabama', 'north dakota', 'minnesota', 'wisconsin', 'massachusetts')
high_states <- c('mississippi', 'arkansas', 'texas', 'alabama')
county_choropleth(mutate(fea, value = POVRATE10), state_zoom = high_low_states, legend = 'Poverty Rate %', title = 'Poverty Rate (2010), Highest and Lowest Food Insecurity States')

#red_color_scale <- c('#ffffb2', '#fed976', '#feb24c', '#fd8d3c', '#fc4e2a', '#e31a1c', '#b10026')
green_color_scale <- c('#edf8fb', '#ccece6', '#99d8c9', '#66c2a4', '#41ae76', '#238b45', '#005824')

county_choropleth(mutate(fea, value = PC_SNAPBEN10), state_zoom = high_low_states, title = 'SNAP benefits per capita (2010), Highest and Lowest Food Insecurity States') + scale_fill_manual(values = green_color_scale, na.value = 'gray50', name = 'Benefits per capita ($)')
## Warning in self$bind(): The following regions were missing and are being
## set to NA: 48033, 48301, 38087, 48269, 25019, 48311, 38007
## Scale for 'fill' is already present. Adding another scale for 'fill',
## which will replace the existing scale.

Let’s zoom into the high food insecurity states

county_choropleth(mutate(fea, value = POVRATE10), state_zoom = high_states, legend = 'Poverty Rate %',
                  title = 'Poverty Rate (2010), High Food Insecurity States')

county_choropleth(mutate(fea, value = PC_SNAPBEN10), state_zoom = high_states, 
                  title = 'SNAP benefits per capita (2010), High Food Insecurity States') + 
  scale_fill_manual(values = green_color_scale, na.value = 'gray50', name = 'Benefits per capita ($)')
## Warning in self$bind(): The following regions were missing and are being
## set to NA: 48033, 48301, 48269, 48311
## Scale for 'fill' is already present. Adding another scale for 'fill',
## which will replace the existing scale.

Let’s look at the state with the worst food insecurity rate, Mississippi.

county_choropleth(mutate(fea, value = POVRATE10), state_zoom = "mississippi", title = 'Poverty Rate (2010) by County for Mississippi')

county_choropleth(mutate(fea, value = PC_SNAPBEN10), state_zoom = "mississippi", title = 'SNAP benefits per capita (2010) for Mississippi')

It looks like SNAP benefits are used in high poverty counties, which means that food assistance money is going to areas that need it.

county_choropleth(mutate(fea, value = SNAPSPTH12), state_zoom = "mississippi", title = 'Number of SNAP-authorized stores/1,000 pop (2012) in Mississippi')

county_choropleth(mutate(fea, value = PCT_18YOUNGER10), state_zoom = "mississippi", 
                  title = '% Population under Age 18 (2010) in Mississippi')

SNAP and WIC Analysis

#df_foodinsec<-read.csv("/home/sanjmeet/food_environ_atlas_analysis/fea_04062017.csv")

# selecting state variables for SNAP
df_SNAP_statevar <- select(fea, FIPS, State, County, PCH_SNAP_09_14, SNAP_PART_RATE08, SNAP_PART_RATE10)

# selecting state variables for WIC
df_WIC_statevar <- select(fea, FIPS, State, County, PCH_WIC_09_14)

# Removing rows(aka counties) which have NA values, since we later want to group by States and be left with single row per state
df_SNAP_statevar <- df_SNAP_statevar[complete.cases(df_SNAP_statevar), ]
df_WIC_statevar <- df_WIC_statevar[complete.cases(df_WIC_statevar), ]

# Finally creating the SNAP and WIC state dataframe. Note, it doesn't matter what we summarize on, Mean/Median/Min/Max; all counties, within the state, have the same value
fea %>% 
  select(-region) %>% 
  group_by(State) %>% 
  summarise_all(mean) -> 
  state_data
## Warning in mean.default(c("02013", "02016", "02020", "02050", "02060",
## "02068", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("01001", "01003", "01005", "01007", "01009",
## "01011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("05001", "05003", "05005", "05007", "05009",
## "05011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("04001", "04003", "04005", "04007", "04009",
## "04011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("06001", "06003", "06005", "06007", "06009",
## "06011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("08001", "08003", "08005", "08007", "08009",
## "08011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("09001", "09003", "09005", "09007", "09009",
## "09011", : argument is not numeric or logical: returning NA
## Warning in mean.default("11001"): argument is not numeric or logical:
## returning NA
## Warning in mean.default(c("10001", "10003", "10005")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("12001", "12003", "12005", "12007", "12009",
## "12011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("13001", "13003", "13005", "13007", "13009",
## "13011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("15001", "15003", "15005", "15007", "15009")):
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("19001", "19003", "19005", "19007", "19009",
## "19011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("16001", "16003", "16005", "16007", "16009",
## "16011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("17001", "17003", "17005", "17007", "17009",
## "17011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("18001", "18003", "18005", "18007", "18009",
## "18011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("20001", "20003", "20005", "20007", "20009",
## "20011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("21001", "21003", "21005", "21007", "21009",
## "21011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("22001", "22003", "22005", "22007", "22009",
## "22011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("25001", "25003", "25005", "25007", "25009",
## "25011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("24001", "24003", "24005", "24009", "24011",
## "24013", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("23001", "23003", "23005", "23007", "23009",
## "23011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("26001", "26003", "26005", "26007", "26009",
## "26011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("27001", "27003", "27005", "27007", "27009",
## "27011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("29001", "29003", "29005", "29007", "29009",
## "29011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("28001", "28003", "28005", "28007", "28009",
## "28011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("30001", "30003", "30005", "30007", "30009",
## "30011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("37001", "37003", "37005", "37007", "37009",
## "37011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("38001", "38003", "38005", "38007", "38009",
## "38011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("31001", "31003", "31005", "31007", "31009",
## "31011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("33001", "33003", "33005", "33007", "33009",
## "33011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("34001", "34003", "34005", "34007", "34009",
## "34011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("35001", "35003", "35005", "35006", "35007",
## "35009", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("32001", "32003", "32005", "32007", "32009",
## "32011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("36001", "36003", "36005", "36007", "36009",
## "36011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("39001", "39003", "39005", "39007", "39009",
## "39011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("40001", "40003", "40005", "40007", "40009",
## "40011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("41001", "41003", "41005", "41007", "41009",
## "41011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("42001", "42003", "42005", "42007", "42009",
## "42011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("44001", "44003", "44005", "44007", "44009")):
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("45001", "45003", "45005", "45007", "45009",
## "45011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("46003", "46005", "46007", "46009", "46011",
## "46013", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("47001", "47003", "47005", "47007", "47009",
## "47011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("48001", "48003", "48005", "48007", "48009",
## "48011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("49001", "49003", "49005", "49007", "49009",
## "49011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("51001", "51003", "51005", "51007", "51009",
## "51011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("50001", "50003", "50005", "50007", "50009",
## "50011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("53001", "53003", "53005", "53007", "53009",
## "53011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("55001", "55003", "55005", "55007", "55009",
## "55011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("54001", "54003", "54005", "54007", "54009",
## "54011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("56001", "56003", "56005", "56007", "56009",
## "56011", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("AK", "AK", "AK", "AK", "AK", "AK", "AK", "AK", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AR", "AR", "AR", "AR", "AR", "AR", "AR", "AR", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AZ", "AZ", "AZ", "AZ", "AZ", "AZ", "AZ", "AZ", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CA", "CA", "CA", "CA", "CA", "CA", "CA", "CA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CO", "CO", "CO", "CO", "CO", "CO", "CO", "CO", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CT", "CT", "CT", "CT", "CT", "CT", "CT", "CT":
## argument is not numeric or logical: returning NA
## Warning in mean.default("DC"): argument is not numeric or logical:
## returning NA
## Warning in mean.default(c("DE", "DE", "DE")): argument is not numeric or
## logical: returning NA
## Warning in mean.default(c("FL", "FL", "FL", "FL", "FL", "FL", "FL", "FL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("GA", "GA", "GA", "GA", "GA", "GA", "GA", "GA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("HI", "HI", "HI", "HI", "HI")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("IA", "IA", "IA", "IA", "IA", "IA", "IA", "IA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ID", "ID", "ID", "ID", "ID", "ID", "ID", "ID", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("IL", "IL", "IL", "IL", "IL", "IL", "IL", "IL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("IN", "IN", "IN", "IN", "IN", "IN", "IN", "IN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("KS", "KS", "KS", "KS", "KS", "KS", "KS", "KS", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("KY", "KY", "KY", "KY", "KY", "KY", "KY", "KY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("LA", "LA", "LA", "LA", "LA", "LA", "LA", "LA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MA", "MA", "MA", "MA", "MA", "MA", "MA", "MA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MD", "MD", "MD", "MD", "MD", "MD", "MD", "MD", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ME", "ME", "ME", "ME", "ME", "ME", "ME", "ME", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MN", "MN", "MN", "MN", "MN", "MN", "MN", "MN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MO", "MO", "MO", "MO", "MO", "MO", "MO", "MO", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MS", "MS", "MS", "MS", "MS", "MS", "MS", "MS", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MT", "MT", "MT", "MT", "MT", "MT", "MT", "MT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ND", "ND", "ND", "ND", "ND", "ND", "ND", "ND", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NE", "NE", "NE", "NE", "NE", "NE", "NE", "NE", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NH", "NH", "NH", "NH", "NH", "NH", "NH", "NH", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NJ", "NJ", "NJ", "NJ", "NJ", "NJ", "NJ", "NJ", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NM", "NM", "NM", "NM", "NM", "NM", "NM", "NM", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NV", "NV", "NV", "NV", "NV", "NV", "NV", "NV", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NY", "NY", "NY", "NY", "NY", "NY", "NY", "NY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OH", "OH", "OH", "OH", "OH", "OH", "OH", "OH", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OK", "OK", "OK", "OK", "OK", "OK", "OK", "OK", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OR", "OR", "OR", "OR", "OR", "OR", "OR", "OR", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("PA", "PA", "PA", "PA", "PA", "PA", "PA", "PA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("RI", "RI", "RI", "RI", "RI")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("SC", "SC", "SC", "SC", "SC", "SC", "SC", "SC", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("SD", "SD", "SD", "SD", "SD", "SD", "SD", "SD", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("TN", "TN", "TN", "TN", "TN", "TN", "TN", "TN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("TX", "TX", "TX", "TX", "TX", "TX", "TX", "TX", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("UT", "UT", "UT", "UT", "UT", "UT", "UT", "UT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("VA", "VA", "VA", "VA", "VA", "VA", "VA", "VA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("VT", "VT", "VT", "VT", "VT", "VT", "VT", "VT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WA", "WA", "WA", "WA", "WA", "WA", "WA", "WA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WI", "WI", "WI", "WI", "WI", "WI", "WI", "WI", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WV", "WV", "WV", "WV", "WV", "WV", "WV", "WV", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WY", "WY", "WY", "WY", "WY", "WY", "WY", "WY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aleutians East", "Aleutians West",
## "Anchorage", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Autauga", "Baldwin", "Barbour", "Bibb",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Arkansas", "Ashley", "Baxter", "Benton",
## "Boone", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Apache", "Cochise", "Coconino", "Gila",
## "Graham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alameda", "Alpine", "Amador", "Butte",
## "Calaveras", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alamosa", "Arapahoe", "Archuleta",
## "Baca", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Fairfield", "Hartford", "Litchfield",
## "Middlesex", : argument is not numeric or logical: returning NA
## Warning in mean.default("District of Columbia"): argument is not numeric or
## logical: returning NA
## Warning in mean.default(c("Kent", "New Castle", "Sussex")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("Alachua", "Baker", "Bay", "Bradford",
## "Brevard", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Appling", "Atkinson", "Bacon", "Baker",
## "Baldwin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Hawaii", "Honolulu", "Kalawao", "Kauai", "Maui":
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Adams", "Allamakee", "Appanoose",
## "Audubon", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Ada", "Adams", "Bannock", "Bear Lake",
## "Benewah", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alexander", "Bond", "Boone", "Brown", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Bartholomew", "Benton",
## "Blackford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allen", "Anderson", "Atchison", "Barber",
## "Barton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Allen", "Anderson", "Ballard",
## "Barren", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Acadia", "Allen", "Ascension", "Assumption", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barnstable", "Berkshire", "Bristol", "Dukes", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allegany", "Anne Arundel", "Baltimore",
## "Calvert", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Androscoggin", "Aroostook", "Cumberland",
## "Franklin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alcona", "Alger", "Allegan", "Alpena",
## "Antrim", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aitkin", "Anoka", "Becker", "Beltrami",
## "Benton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Andrew", "Atchison", "Audrain",
## "Barry", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alcorn", "Amite", "Attala", "Benton", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaverhead", "Big Horn", "Blaine",
## "Broadwater", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alamance", "Alexander", "Alleghany", "Anson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Barnes", "Benson", "Billings",
## "Bottineau", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Antelope", "Arthur", "Banner",
## "Blaine", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Belknap", "Carroll", "Cheshire", "Coos",
## "Grafton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Atlantic", "Bergen", "Burlington", "Camden", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bernalillo", "Catron", "Chaves", "Cibola",
## "Colfax", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Churchill", "Clark", "Douglas", "Elko",
## "Esmeralda", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Allegany", "Bronx", "Broome",
## "Cattaraugus", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Ashland", "Ashtabula",
## "Athens", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Alfalfa", "Atoka", "Beaver",
## "Beckham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Baker", "Benton", "Clackamas", "Clatsop",
## "Columbia", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allegheny", "Armstrong", "Beaver",
## "Bedford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bristol", "Kent", "Newport", "Providence",
## "Washington": argument is not numeric or logical: returning NA
## Warning in mean.default(c("Abbeville", "Aiken", "Allendale", "Anderson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aurora", "Beadle", "Bennett", "Bon Homme",
## "Brookings", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Bedford", "Benton", "Bledsoe",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Andrews", "Angelina", "Aransas", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaver", "Box Elder", "Cache", "Carbon",
## "Daggett", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Accomack", "Albemarle", "Alleghany", "Amelia", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Addison", "Bennington", "Caledonia",
## "Chittenden", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Asotin", "Benton", "Chelan",
## "Clallam", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Ashland", "Barron", "Bayfield",
## "Brown", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barbour", "Berkeley", "Boone", "Braxton",
## "Brooke", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Big Horn", "Campbell", "Carbon",
## "Converse", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("AK", "AK", "AK", "AK", "AK", "AK", "AK", "AK", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AR", "AR", "AR", "AR", "AR", "AR", "AR", "AR", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("AZ", "AZ", "AZ", "AZ", "AZ", "AZ", "AZ", "AZ", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CA", "CA", "CA", "CA", "CA", "CA", "CA", "CA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CO", "CO", "CO", "CO", "CO", "CO", "CO", "CO", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("CT", "CT", "CT", "CT", "CT", "CT", "CT", "CT":
## argument is not numeric or logical: returning NA
## Warning in mean.default("DC"): argument is not numeric or logical:
## returning NA
## Warning in mean.default(c("DE", "DE", "DE")): argument is not numeric or
## logical: returning NA
## Warning in mean.default(c("FL", "FL", "FL", "FL", "FL", "FL", "FL", "FL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("GA", "GA", "GA", "GA", "GA", "GA", "GA", "GA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("HI", "HI", "HI", "HI", "HI")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("IA", "IA", "IA", "IA", "IA", "IA", "IA", "IA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ID", "ID", "ID", "ID", "ID", "ID", "ID", "ID", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("IL", "IL", "IL", "IL", "IL", "IL", "IL", "IL", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("IN", "IN", "IN", "IN", "IN", "IN", "IN", "IN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("KS", "KS", "KS", "KS", "KS", "KS", "KS", "KS", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("KY", "KY", "KY", "KY", "KY", "KY", "KY", "KY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("LA", "LA", "LA", "LA", "LA", "LA", "LA", "LA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MA", "MA", "MA", "MA", "MA", "MA", "MA", "MA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MD", "MD", "MD", "MD", "MD", "MD", "MD", "MD", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ME", "ME", "ME", "ME", "ME", "ME", "ME", "ME", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MN", "MN", "MN", "MN", "MN", "MN", "MN", "MN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MO", "MO", "MO", "MO", "MO", "MO", "MO", "MO", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MS", "MS", "MS", "MS", "MS", "MS", "MS", "MS", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("MT", "MT", "MT", "MT", "MT", "MT", "MT", "MT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("ND", "ND", "ND", "ND", "ND", "ND", "ND", "ND", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NE", "NE", "NE", "NE", "NE", "NE", "NE", "NE", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NH", "NH", "NH", "NH", "NH", "NH", "NH", "NH", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NJ", "NJ", "NJ", "NJ", "NJ", "NJ", "NJ", "NJ", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NM", "NM", "NM", "NM", "NM", "NM", "NM", "NM", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NV", "NV", "NV", "NV", "NV", "NV", "NV", "NV", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("NY", "NY", "NY", "NY", "NY", "NY", "NY", "NY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OH", "OH", "OH", "OH", "OH", "OH", "OH", "OH", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OK", "OK", "OK", "OK", "OK", "OK", "OK", "OK", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("OR", "OR", "OR", "OR", "OR", "OR", "OR", "OR", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("PA", "PA", "PA", "PA", "PA", "PA", "PA", "PA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("RI", "RI", "RI", "RI", "RI")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("SC", "SC", "SC", "SC", "SC", "SC", "SC", "SC", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("SD", "SD", "SD", "SD", "SD", "SD", "SD", "SD", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("TN", "TN", "TN", "TN", "TN", "TN", "TN", "TN", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("TX", "TX", "TX", "TX", "TX", "TX", "TX", "TX", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("UT", "UT", "UT", "UT", "UT", "UT", "UT", "UT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("VA", "VA", "VA", "VA", "VA", "VA", "VA", "VA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("VT", "VT", "VT", "VT", "VT", "VT", "VT", "VT", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WA", "WA", "WA", "WA", "WA", "WA", "WA", "WA", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WI", "WI", "WI", "WI", "WI", "WI", "WI", "WI", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WV", "WV", "WV", "WV", "WV", "WV", "WV", "WV", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("WY", "WY", "WY", "WY", "WY", "WY", "WY", "WY", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aleutians East", "Aleutians West",
## "Anchorage", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Autauga", "Baldwin", "Barbour", "Bibb",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Arkansas", "Ashley", "Baxter", "Benton",
## "Boone", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Apache", "Cochise", "Coconino", "Gila",
## "Graham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alameda", "Alpine", "Amador", "Butte",
## "Calaveras", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alamosa", "Arapahoe", "Archuleta",
## "Baca", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Fairfield", "Hartford", "Litchfield",
## "Middlesex", : argument is not numeric or logical: returning NA
## Warning in mean.default("District of Columbia"): argument is not numeric or
## logical: returning NA
## Warning in mean.default(c("Kent", "New Castle", "Sussex")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("Alachua", "Baker", "Bay", "Bradford",
## "Brevard", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Appling", "Atkinson", "Bacon", "Baker",
## "Baldwin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Hawaii", "Honolulu", "Kalawao", "Kauai", "Maui":
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Adams", "Allamakee", "Appanoose",
## "Audubon", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Ada", "Adams", "Bannock", "Bear Lake",
## "Benewah", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alexander", "Bond", "Boone", "Brown", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Bartholomew", "Benton",
## "Blackford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allen", "Anderson", "Atchison", "Barber",
## "Barton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Allen", "Anderson", "Ballard",
## "Barren", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Acadia", "Allen", "Ascension", "Assumption", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barnstable", "Berkshire", "Bristol", "Dukes", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allegany", "Anne Arundel", "Baltimore",
## "Calvert", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Androscoggin", "Aroostook", "Cumberland",
## "Franklin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alcona", "Alger", "Allegan", "Alpena",
## "Antrim", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aitkin", "Anoka", "Becker", "Beltrami",
## "Benton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Andrew", "Atchison", "Audrain",
## "Barry", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alcorn", "Amite", "Attala", "Benton", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaverhead", "Big Horn", "Blaine",
## "Broadwater", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alamance", "Alexander", "Alleghany", "Anson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Barnes", "Benson", "Billings",
## "Bottineau", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Antelope", "Arthur", "Banner",
## "Blaine", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Belknap", "Carroll", "Cheshire", "Coos",
## "Grafton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Atlantic", "Bergen", "Burlington", "Camden", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bernalillo", "Catron", "Chaves", "Cibola",
## "Colfax", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Churchill", "Clark", "Douglas", "Elko",
## "Esmeralda", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Allegany", "Bronx", "Broome",
## "Cattaraugus", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Ashland", "Ashtabula",
## "Athens", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Alfalfa", "Atoka", "Beaver",
## "Beckham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Baker", "Benton", "Clackamas", "Clatsop",
## "Columbia", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allegheny", "Armstrong", "Beaver",
## "Bedford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bristol", "Kent", "Newport", "Providence",
## "Washington": argument is not numeric or logical: returning NA
## Warning in mean.default(c("Abbeville", "Aiken", "Allendale", "Anderson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aurora", "Beadle", "Bennett", "Bon Homme",
## "Brookings", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Bedford", "Benton", "Bledsoe",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Andrews", "Angelina", "Aransas", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaver", "Box Elder", "Cache", "Carbon",
## "Daggett", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Accomack", "Albemarle", "Alleghany", "Amelia", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Addison", "Bennington", "Caledonia",
## "Chittenden", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Asotin", "Benton", "Chelan",
## "Clallam", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Ashland", "Barron", "Bayfield",
## "Brown", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barbour", "Berkeley", "Boone", "Braxton",
## "Brooke", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Big Horn", "Campbell", "Carbon",
## "Converse", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aleutians East", "Aleutians West",
## "Anchorage", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Autauga", "Baldwin", "Barbour", "Bibb",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Arkansas", "Ashley", "Baxter", "Benton",
## "Boone", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Apache", "Cochise", "Coconino", "Gila",
## "Graham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alameda", "Alpine", "Amador", "Butte",
## "Calaveras", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alamosa", "Arapahoe", "Archuleta",
## "Baca", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Fairfield", "Hartford", "Litchfield",
## "Middlesex", : argument is not numeric or logical: returning NA
## Warning in mean.default("District of Columbia"): argument is not numeric or
## logical: returning NA
## Warning in mean.default(c("Kent", "New Castle", "Sussex")): argument is not
## numeric or logical: returning NA
## Warning in mean.default(c("Alachua", "Baker", "Bay", "Bradford",
## "Brevard", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Appling", "Atkinson", "Bacon", "Baker",
## "Baldwin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Hawaii", "Honolulu", "Kalawao", "Kauai", "Maui":
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Adams", "Allamakee", "Appanoose",
## "Audubon", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Ada", "Adams", "Bannock", "Bear Lake",
## "Benewah", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alexander", "Bond", "Boone", "Brown", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Bartholomew", "Benton",
## "Blackford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allen", "Anderson", "Atchison", "Barber",
## "Barton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Allen", "Anderson", "Ballard",
## "Barren", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Acadia", "Allen", "Ascension", "Assumption", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barnstable", "Berkshire", "Bristol", "Dukes", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Allegany", "Anne Arundel", "Baltimore",
## "Calvert", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Androscoggin", "Aroostook", "Cumberland",
## "Franklin", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alcona", "Alger", "Allegan", "Alpena",
## "Antrim", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aitkin", "Anoka", "Becker", "Beltrami",
## "Benton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Andrew", "Atchison", "Audrain",
## "Barry", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Alcorn", "Amite", "Attala", "Benton", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaverhead", "Big Horn", "Blaine",
## "Broadwater", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Alamance", "Alexander", "Alleghany", "Anson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Barnes", "Benson", "Billings",
## "Bottineau", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Antelope", "Arthur", "Banner",
## "Blaine", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Belknap", "Carroll", "Cheshire", "Coos",
## "Grafton", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Atlantic", "Bergen", "Burlington", "Camden", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bernalillo", "Catron", "Chaves", "Cibola",
## "Colfax", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Churchill", "Clark", "Douglas", "Elko",
## "Esmeralda", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Allegany", "Bronx", "Broome",
## "Cattaraugus", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allen", "Ashland", "Ashtabula",
## "Athens", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adair", "Alfalfa", "Atoka", "Beaver",
## "Beckham", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Baker", "Benton", "Clackamas", "Clatsop",
## "Columbia", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Allegheny", "Armstrong", "Beaver",
## "Bedford", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Bristol", "Kent", "Newport", "Providence",
## "Washington": argument is not numeric or logical: returning NA
## Warning in mean.default(c("Abbeville", "Aiken", "Allendale", "Anderson", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Aurora", "Beadle", "Bennett", "Bon Homme",
## "Brookings", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Bedford", "Benton", "Bledsoe",
## "Blount", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Anderson", "Andrews", "Angelina", "Aransas", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Beaver", "Box Elder", "Cache", "Carbon",
## "Daggett", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Accomack", "Albemarle", "Alleghany", "Amelia", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("Addison", "Bennington", "Caledonia",
## "Chittenden", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Asotin", "Benton", "Chelan",
## "Clallam", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Adams", "Ashland", "Barron", "Bayfield",
## "Brown", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Barbour", "Berkeley", "Boone", "Braxton",
## "Brooke", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("Albany", "Big Horn", "Campbell", "Carbon",
## "Converse", : argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default("1.000000"): argument is not numeric or logical:
## returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000")): argument is
## not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "1.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "0.000000", "1.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "1.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "1.000000", "1.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "1.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "1.000000", "1.000000", "1.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "1.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "1.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("1.000000", "0.000000", "1.000000", "0.000000", :
## argument is not numeric or logical: returning NA
## Warning in mean.default(c("0.000000", "0.000000", "0.000000", "1.000000", :
## argument is not numeric or logical: returning NA
df_SNAP_statevar %>% 
  group_by(State) %>% 
  summarise_all(min) -> 
  df_SNAP_state

 df_WIC_statevar %>% 
  group_by(State) %>% 
  summarise_all(min) -> 
   df_WIC_state

g1 <- ggplot(df_SNAP_state, aes(x = State, y = df_SNAP_state$PCH_SNAP_09_14 / 100)) + 
  geom_bar(stat = "identity", color = "black", fill = "lightblue") + 
  scale_y_continuous(labels = scales::percent) + 
  coord_flip() + 
  ylab("Change") + 
  xlab("State") + 
  ggtitle("Percent change in SNAP participation from 2009 to 2014")

g2 <- ggplot(df_WIC_state, aes(x = State, y = df_WIC_state$PCH_WIC_09_14 / 100)) + 
  geom_bar(stat = "identity", color = "black", fill = "lightpink") + 
  scale_y_continuous(labels = scales::percent) + 
  coord_flip() + 
  ylab("Change") + 
  xlab("State") + 
  ggtitle("Percent change in WIC participation from 2009 to 2014")

grid.arrange(g1, g2, ncol = 2)

##### Dodge columns: plot both SNAP and WIC per state, and sort by SNAP number

As mentioned by numerous reports, both by the USDA and the general media, SNAP popularity/participation increased from 2009-2014. Nevada was the only state for which this wasn’t true.

On the other hand, WIC participation decreased across all states, except in New York and Minnesota.

This can mainly be attributed to the following reasons:

  1. SNAP, unlike WIC, provides dollar amount that is easier to spend and manage than coupons.

  2. SNAP registration process is easier and does not require a certificate of approval from a physician. WIC registration,on the other hand, is longer and more rigorous.

  3. SNAP applies to all members of a family, while WIC applies to Women, Infants and Children.

  4. A decrease in the overall birth rate across all states.

# Lets dive deeper into the previous findings by taking a closer look at the SNAP and WIC redemptions in authorized stores in 2008 and 2012. We will look at the data of five states with high food insecurity and five states with low food insecurity here.

# Four Target States with highest food insecurity[according to USDA 2013-2014]: Mississippi,Arkansas, Texas, Alabama
# Four Target states with lowest food insecurity[according to USDA 2013-2014]: North Dakota,Minnesota, Wisconsin, Massachusetts

#selecting county level variables for SNAP
df_SNAP_countyvar <- select(fea, FIPS, State, County, SNAPSPTH08, SNAPSPTH12, PCH_SNAPSPTH_08_12, REDEMP_SNAPS08, REDEMP_SNAPS12, PCH_REDEMP_SNAPS_08_12, PC_SNAPBEN08, PC_SNAPBEN10, PCH_PC_SNAPBEN_08_10)

#Selecting county level variables for WIC
df_WIC_countyvar <- select(fea, FIPS, State, County, PC_WIC_REDEMP08, PC_WIC_REDEMP12, REDEMP_WICS08, REDEMP_WICS12)

#removing counties for which SNAP information is not available
df_SNAP_countyvar <- df_SNAP_countyvar[complete.cases(df_SNAP_countyvar), ]

#removing counties for which WIC information is not available
df_WIC_countyvar <- df_WIC_countyvar[complete.cases(df_WIC_countyvar), ]

# keeping counties for only target states
df_SNAP_targetcounty <- filter(df_SNAP_countyvar, State %in% c(high_insecurity_states, low_insecurity_states))
df_WIC_targetcounty <- filter(df_WIC_countyvar, State %in% c(high_insecurity_states, low_insecurity_states))

#Adding new columns to contain redemption values in 1k scale
df_SNAP_targetcounty <- df_SNAP_targetcounty %>% 
  mutate(Redemp08_1K = REDEMP_SNAPS08 / 1000.00, Redemp12_1K = REDEMP_SNAPS12 / 1000.00)
df_WIC_targetcounty <- df_WIC_targetcounty %>% 
  mutate(Redemp08_1K = REDEMP_WICS08 / 1000.00, Redemp12_1K = REDEMP_WICS12 / 1000.00)

g1 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +Redemp08_1K), y = df_SNAP_targetcounty$Redemp08_1K)) + 
  geom_boxplot(fill = "azure") + 
  xlab("State") + 
  ylab("Average Redemption (in $1k)") + 
  ggtitle("Average Monthly SNAP Redemptions in Authorized Stores in 2008") + 
  ylim(0, 600) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$Redemp08_1K, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 4.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
g3 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +Redemp12_1K), y = df_SNAP_targetcounty$Redemp12_1K)) + 
  geom_boxplot() + 
  xlab("State") + 
  ylab("Average redemption(in $1k)") + 
  ggtitle("Average Monthly SNAP Redemptions in Authorized Stores in 2012") + 
  ylim(0, 600) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$Redemp12_1K, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 2.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 2.5, xmax = 3.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  annotate('rect', xmin = 3.5, xmax = 5.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 5.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
g2 <- ggplot(df_WIC_targetcounty, aes(x = reorder(State, +Redemp08_1K), y = df_WIC_targetcounty$Redemp08_1K)) + 
  geom_boxplot(fill = "#c2a5cf") + 
  xlab("State") + 
  ylab("Average Redemption (in $1k)") + 
  ggtitle("Average Monthly WIC Redemptions in Authorized Stores in 2008") + 
  ylim(0, 600) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_WIC_targetcounty$Redemp08_1K, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 4.5, xmax = 7.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
g4 <- ggplot(df_WIC_targetcounty, aes(x = reorder(State, +Redemp12_1K), y = df_WIC_targetcounty$Redemp12_1K)) + 
  geom_boxplot() + 
  xlab("State") + 
  ylab("Average Redemption (in $1k)") + 
  ggtitle("Average Monthly WIC Redemptions in Authorized Stores in 2012") + 
  ylim(0, 600) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_WIC_targetcounty$Redemp12_1K, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 4.5, xmax = 7.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  geom_boxplot(fill = "ghostwhite") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
grid.arrange(g1, g2, g3, g4, ncol = 2)

##### Make sure all y-axes are 0 - 600
##### Does the comparison make sense?

Key observations:

  1. Average monthly SNAP redemptions increased in across all states from 2008 to 2012. The US median also increased.
  2. No clear pattern for WIC redemptions. The Us median stayed the same. The states either see a decline or remain the same. Unlike SNAP, none of the states see an increase.
  3. Both in 2008 and 2012, the US median for average monthly SNAP redemptions was significantly higher than the average monthly WIC redemptions.
  4. Also, as should be the case, the average redemptions were higher in states with higher food insecurity.

Since SNAP seems to be the popular choice of food assistance program, we decided to explore some of the other state and county variables associated with it.

We earlier saw how SNAP participation increased, across all US states (except Nevada) from 2008 to 2012. We now look at SNAP participation as a percentage of eligible participants in the year 2008 and 2012. This shall give us a more concrete proof of whether the participation actually increased or not.

# 2008
g1 <- ggplot(df_SNAP_state, aes(x = State, y = df_SNAP_state$SNAP_PART_RATE08 / 100)) + 
  geom_point() + 
  scale_y_continuous(labels = scales::percent) + 
  coord_flip() + 
  ylab("SNAP participants (% eligible pop)") + 
  xlab("State") + 
  ggtitle("SNAP participation as a % of eligible participants in 2008, with US median") + 
  geom_hline(yintercept = median(df_SNAP_state$SNAP_PART_RATE08 / 100, na.rm = TRUE), linetype = 'longdash')

# 2010
g2 <- ggplot(df_SNAP_state, aes(x = State, y = df_SNAP_state$SNAP_PART_RATE10 / 100)) + 
  geom_point() + 
  scale_y_continuous(labels = scales::percent) + 
  coord_flip() + 
  ylab("SNAP participants (% eligible pop)") + 
  xlab("State") + 
  ggtitle("SNAP participation as a % of eligible population in 2010, with US median") + 
  geom_hline(yintercept = median(df_SNAP_state$SNAP_PART_RATE10 / 100, na.rm = TRUE), linetype = 'longdash')

grid.arrange(g1, g2, ncol = 1)

##### Use Cleveland dot plot style from above, coloring by year; 2 lines (labeled, with color); reorder

Key observations:

  1. The median SNAP participation as a percentage of eligible populations increased from 2008 to 2012.
  2. Even for Nevada, the participation as a percentage of eligible population increased from 2008 to 2012.

Next, we look at SNAP authorized stores per 1000 people at county level for our five food secure and five food insecure states.

# County Level Exploration

# Boxplot for SNAP authorized stores/1000 population in 2008 across counties of these 9 states
g1 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +SNAPSPTH08), y = SNAPSPTH08)) + 
  geom_boxplot(fill = "#a6dba0") + 
  xlab("State") + 
  ylab("Stores to Population Ratio") + 
  ggtitle("SNAP Authorized Stores per 1000 people in 2008") + 
  ylim(0, 3) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$SNAPSPTH08, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 3.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 3.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  annotate('rect', xmin = 4.5, xmax = 5.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 5.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()

#Boxplot for SNAP authorised stores/1000 population in 2012 across counties of these 9 states
g2 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +SNAPSPTH12), y = SNAPSPTH12)) + 
  geom_boxplot(fill = "#dcf1da") + 
  xlab("State") + 
  ggtitle("SNAP Authorized Stores per 1000 people in 2012") + 
  ylim(0,3) + 
  theme(axis.title.y = element_blank()) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$SNAPSPTH12, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 3.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 3.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  annotate('rect', xmin = 4.5, xmax = 5.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 5.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()

grid.arrange(g1, g2, ncol = 2)

Observe how the number of stores/1000 people is higher in states with higher food insecurity.

Notice how the number of stores/1000 population increased nearly everywhere and, consequently, the US median. This is something that would have been expected given recession and incoming of the Obama government.

For the recession period from 2008-2013, the US government had announced that the SNAP benefit amount had been increased , in all states, to help combat food insecurity due to unemployment. Let’s see if the data can confirm this fact. We look at the SNAP Benefits per capita variable for our chosen 10 states.

# Lastly we look at SNAP benefits per capita

g1 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +PC_SNAPBEN08), y = PC_SNAPBEN08)) + 
  geom_boxplot(fill = "lightblue") + 
  xlab("State") + 
  ylab("Avg benefits per capita") + 
  ggtitle("Average Monthly SNAP Benefits per capita in 2008") + 
  ylim(0, 60) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$PC_SNAPBEN08, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 4.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
g2 <- ggplot(df_SNAP_targetcounty, aes(x = reorder(State, +PC_SNAPBEN10), y = PC_SNAPBEN10)) + 
  geom_boxplot(fill = "lightblue") + 
  xlab("State") + 
  ylab("Average monthly SNAP benefits per capita") + 
  ggtitle("Average Monthly SNAP Benefits per capita in 2010") + 
  ylim(0, 60) + 
  scale_y_continuous(labels = scales::dollar) + 
  geom_hline(yintercept = median(df_SNAP_targetcounty$PC_SNAPBEN10, na.rm = TRUE), linetype = 'longdash') + 
  annotate('rect', xmin = 0.5, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.2, fill = "cyan3") + 
  annotate('rect', xmin = 4.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.1, fill = "red") + 
  theme_bw()
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.
grid.arrange(g1, g2, ncol = 2)

It is true, the SNAP benefit amount was increased across all states from 2008-2012.

The states with higher food insecurity seem to have much a higher SNAP benefits per capita amount vs the states with lower food insecurity.

In the previous graph we saw how SNAP benefits per capita increased across all states from 2008 to 2012. Although this was reassuring, we wanted to see at how uniform this increase was across all the counties of a highly food insecure state. We chose to look at the counties of Mississipi.

df_benefits_perchange <- select(df_SNAP_targetcounty, FIPS, State, County, PCH_PC_SNAPBEN_08_10)
df_MSbenperchange <- filter(df_benefits_perchange, State == "MS")
df_MSbenperchange <- df_MSbenperchange %>% mutate(Year = 2010)
colnames(df_MSbenperchange)[4] <- "Percentage_Change"

ggplot(df_MSbenperchange, aes(Year, County, fill = Percentage_Change)) + 
  geom_tile() + 
  scale_fill_viridis() + 
  ggtitle("Percent Change in SNAP Benefits per capita in the counties of Mississipi") + 
  ylab("Country") + 
  xlab("2008-2010") + 
  theme(axis.text.x = element_blank(), axis.ticks.x = element_blank())

As can be seen from the Heat Map above, the increase in SNAP benefit amount was not at all uniform. There were significant number of states which saw less than a 60% increase in the benefit amount. Mississipi being the highest food insecure state, we had expected to see more counties with a higher increase in the SNAP benefit amount.